
/*
 * 1-18-92      Jack Snoeyink
 * Declarations for simple geometric algorithms
 * (under construction)
 */

#ifndef GEOM_DEFS

#define GEOM_DEFS

#include <math.h>

#ifndef FALSE
#define FALSE 0
#endif
#ifndef TRUE
#define TRUE 1
#endif
#ifndef NULL
#define NULL 0
#endif

#define TWO_PI 6.283185
#define COORDTYPE double
#define GEOMETRY HOMOG
#define VECTORSIZE 4

#define vH(p) v4d((COORDTYPE *) &(p)); /* for SGI graphics */
#define vHp(p) v4d((COORDTYPE *) (p));

typedef struct _POINT {		/* compatible with SGI homog coords (yech)  */
  COORDTYPE x, y, z, w;
  COORDTYPE el;			/* elevation */
} POINT;

#define RandomPOINT(p, t) /* Random point uniform over [-t, t]^3 */\
 {(p)->w = 1.0;\
  (p)->x = t*(2.0 * drand48() - 1);\
  (p)->y = t*(2.0 * drand48() - 1);\
  (p)->z = t*(2.0 * drand48() - 1);}

#define LINCOMB(a, p, b, q,  r) /* 3-d homog vector linear combination */\
  (r)->w = (p)->w*(q)->w; \
  (r)->x = (a)*(p)->x*(q)->w + (b)*(q)->x*(p)->w; \
  (r)->y = (a)*(p)->y*(q)->w + (b)*(q)->y*(p)->w; \
  (r)->z = (a)*(p)->z*(q)->w + (b)*(q)->z*(p)->w; 

#define LINCOMB_C(a, p, b, q,  r) /* Cartesian lincomb (assumes p.w=q.w) */\
  (r)->w = (p)->w; \
  (r)->x = (a)*(p)->x + (b)*(q)->x; \
  (r)->y = (a)*(p)->y + (b)*(q)->y; \
  (r)->z = (a)*(p)->z + (b)*(q)->z; 

#define OPERATE_C(r, op, p)	/* 3-d cartesian vector unary operation */\
  (r)->w = (p)->w; \
  (r)->x op (p)->x; \
  (r)->y op (p)->y; \
  (r)->z op (p)->z;


#define DOT(p,q)		/* 3-d homog dot prod */\
  ((p)->w*(q)->w + (p)->x*(q)->x + (p)->y*(q)->y + (p)->z*(q)->z)

/* Determinants by minor expansion */

#define DET2x2(p, q, i, j) \
  ((p)->i*(q)->j - (p)->j*(q)->i)

#define DET3x3(p, q, r, i, j, k) \
  (  (p)->i*DET2x2(q,r,j,k) \
   - (p)->j*DET2x2(q,r,i,k) \
   + (p)->k*DET2x2(q,r,i,j))

#define MEET(p,q,r, s)		/* 3-d homog vector meet */\
  (s)->w =   DET3x3(p, q, r, x, y, z);\
  (s)->x = - DET3x3(p, q, r, w, y, z);\
  (s)->y =   DET3x3(p, q, r, w, x, z);\
  (s)->z = - DET3x3(p, q, r, w, x, y);

#ifndef EPSILON  
#define EPSILON 1.0e-9
#endif
#define EPSILONSQ 1.0e-16

#define NegativeP(x) (x < -EPSILON)
#define PositiveP(x) (x > EPSILON)
#define ZeroP(x) (fabs(x) <= EPSILON) /* single evaluation */
#define FastZeroP(x) ((-EPSILON <= x) && (x <= EPSILON))

#define SGN(x) (NegativeP(x) ? -1 : (PositiveP(x) ? 1 : 0))

#define ON_PL(p, l) ZeroP(DOT(p, l))
#define OFF_PL(p, l) !ZeroP(DOT(p, l))
#define UP_PL(p, l) PositiveP(DOT(p, l))
#define DOWN_PL(p, l) NegativeP(DOT(p, l))
#define UP_ON_PL(p, l) !NegativeP(DOT(p, l))
#define DOWN_ON_PL(p, l) !PositiveP(DOT(p, l))

#endif